[1]:


from IPython.display import IFrame

# Youtube
IFrame("https://www.youtube.com/embed/D8S4GGmHUMw", "100%", 500)


[1]:

Configurações

[2]:


import numpy as np
from scipy.integrate import odeint

from bokeh.palettes import brewer, Inferno10
from bokeh.plotting import figure, show
from bokeh.io import output_notebook
from bokeh.layouts import gridplot

output_notebook()
Loading BokehJS ...

Modelo para pequeno grupo

Há no senso comum o entendimento que há uma chance de “pegar” uma doença. Supondo que a probailidade que representa esta sensação seja \(\mathbf{p}\) e que a probilidade de escapar desta doença seja \(\mathbf{q}\). Vamos admitir um cenário onde tenhamos uma família de três pessoas \(\mathbf{N}=3\) e um dos individuos estaje infectado (I). Supondo que neste caso \(\mathbf{p}=0.4\) logo \(\mathbf{q}=0.6\). Admitindo este cenário pode-se explorar as seguintes progressões da doença (\(P_i\) - infectado e \(\overline{P_i}\) - não infectado) :

Pensando em casos:

  1. Apenas o primeiro ifectado desenvovle a doença.
  2. O infectado transmite para um dos outros dois suscetíveis apenas.
  3. Os dois sucetíveis são infectados um após o outro.
  4. Ambos os suscetíveis são infectados simultaneamente.

Traduzindo em cadeias de transmissão:

  1. 1 \(\rightarrow\) 0
  2. 1 \(\rightarrow\) 1 \(\rightarrow\) 0
  3. 1 \(\rightarrow\) 1 \(\rightarrow\) 1
  4. 1 \(\rightarrow\) 2

Traduzindo em probabilidades de cada cadeia de transmissão:

  1. \(\overline{P_1}\) e \(\overline{P_2}\) - \(\mathbf{q^2}\)
  2. \(\overline{P_1}\) e \(P_2\) - \(\mathbf{pq^2}\) ou \(P_1\) e \(\overline{P_2}\) - \(\mathbf{pq^2}\)
  3. \(\overline{P_1}\) e \(P_2\) \(\rightarrow\) \(P_1\) e \(P_2\) - \(\mathbf{p^2q}\) ou \(P_1\) e \(\overline{P_2}\) \(\rightarrow\) \(P_1\) e \(P_2\) - \(\mathbf{p^2q}\)
  4. \(P_1\) e \(P_2\) - \(\mathbf{p^2}\)
[3]:


p = 0.4
q = 1 - p

C1 = q**2
C2 = 2*p*(q**2)
C3 = 2*q*(p**2)
C4 = p**2

print(C1,C2,C3,C4, C3+C4, C1+C2+C3+C4)


0.36 0.288 0.19200000000000003 0.16000000000000003 0.3520000000000001 1.0

Considerando as probabilidades calculadas para cada uma das sequências, aplicada um cenário de 1000 famílias com três integrantes, pode-se estimar que cerca de 360 famílias, atigidas pela doença, terão apenas 1 infectado, 288 famílias terão pelo menos 1 infectado e 352 famílias todos serão infectados.

Se generalizarmos a análise pode-se dizer que o número de infectados \(\mathbf{I_{t}}\) na próxima geração de infectados será \(\mathbf{I_{t+1}}\). O termo \(\mathbf{t}\) indica a geração da epidemia e indica com uma certa imprecisão a escala de tempo.

\begin{align*} P(I_{t+1}=i_{t+1}\vert S_t=s_t,I_t=i_t) \sim \binom{s_t}{i_t} (1-q^i_t)^{I_{t+1}} (q)^{i_t(s_t-i_{t+1})}, s_t \geq i_{t+1} \end{align*}

Generalizando o modelo

A partir da proposta que o contágio é dependente de uma probabilidade \(\mathbf{p}\) de transmissão da doença, que corresponde a uma probabilidade \(\mathbf{q}\) de escapar de ser contaminado, pode-se estabelecer um modelo que seja capaz de representar a cadeia de transmissão.

\[\begin{split}\begin{split} I_{t+1} & \sim binomial(S(t),p) \\ S_{t+1} & = S_t - I_{t+1} \\ R_{t+1} & = R_t + I_t \end{split}\end{split}\]

Observe a semelhança com o modelo determinístico, exceto o fato do número de infectados ser um valor randômico, governado pela distribuição Binomial.

[4]:


# Modelo Reed-Frost
ngen  = 30;     # número de gerações
Sinit = 2000;   # população suscetível
Iinit = 1;      # infectados
Rinit = 0;      # removidos / recuperados

q = 0.999;

nsims = 2000;   # número de simulações

x, y = list(), list()
M = np.zeros((nsims-1,ngen-1))

TOOLS="zoom_in,zoom_out,save"
p = figure(tools=TOOLS, plot_width=600, plot_height=400)

for i in range(1,nsims-1):
    S = np.linspace(0,1,ngen) - np.linspace(0,1,ngen)
    I = np.linspace(0,1,ngen) - np.linspace(0,1,ngen)
    R = np.linspace(0,1,ngen) - np.linspace(0,1,ngen)
    S[1], I[1], R[1] = Sinit, Iinit, Rinit
    for j in range(1,ngen-1):
        # np.random.binomial(n, p, 1000)
        I[j+1] = np.random.binomial(S[j], 1 - q**I[j], 1);
        S[j+1] = S[j] - I[j+1]
        R[j+1] = R[j] + I[j]
        M[i,j] = I[j]

    p.line(range(ngen), I, line_width=2, color="#8e44ad", line_alpha=0.05)

p.grid.grid_line_alpha = 0
p.ygrid.band_fill_color = "olive"
p.ygrid.band_fill_alpha = 0.1
p.yaxis.axis_label = "Infectados"
p.xaxis.axis_label = "Gerações"

show(p)

Análise dos cenários para cada geração nova de infectados

Observe que para cada nova geração a distribuição estatística varia e pode-se observar o comportamento da distribuição Binomial. Neste sentido pode-se observar quem “governa” a incerteza da evolução da epidemia. A linha pontilhada em cada histograma representa a média para cada geração.

[5]:


num = 10
x = np.linspace(0, 100, 1000)
pallete = Inferno10

fig_list = []
for i in range(num):
    hist, edges = np.histogram(M[:,i], density=False, bins=50)
    p = figure(tools="pan,hover,lasso_select", plot_width=300, plot_height=300)
    p.quad(top=hist, bottom=0, left=edges[:-1], right=edges[1:],
     fill_color=pallete[i], fill_alpha=0.5, line_color=pallete[i],
     legend_label="k = " + str(i) )
    p.grid.grid_line_alpha = 0
    p.ygrid.band_fill_alpha = 0.1
    p.ygrid.band_fill_color = "olive"
    p.yaxis.axis_label = "Frequência"
    p.xaxis.axis_label = "Infectados"
    fig_list.append(p)

grid = gridplot(fig_list, ncols=2)
show(grid)


[6]:


num = 15

p = figure(tools="pan,hover,lasso_select", y_axis_type="log", plot_width=600, plot_height=400)

for i in range(num):
    x = (i+1) * np.ones((len(M[:,i]),))
    try:
        p.scatter(x, M[:,i], size=5, fill_color=pallete[i], fill_alpha=0.05, line_alpha=0)
    except:
        p.scatter(x, M[:,i], size=5, fill_color=pallete[i-10], fill_alpha=0.05, line_alpha=0)
p.grid.grid_line_alpha = 0
p.ygrid.band_fill_alpha = 0.1
p.ygrid.band_fill_color = "olive"
p.yaxis.axis_label = "Infectados"
p.xaxis.axis_label = "Gerações"

show(p)

Evolução da média e desvio padrão

Para uma visão mais clara podemos visualizar a média e o desvio padrão para cada geração.

[7]:

a = M.mean(0)
e = M.std(0)
x = np.linspace(0,ngen, ngen-1)

p = figure(tools="pan,hover,lasso_select", plot_width=600, plot_height=400)

p.varea(x=x, y1=a+e/2, y2=a-e/2, legend_label="Desvio Padrão", fill_alpha=0.3, fill_color="#3498db")
p.line(x, a, legend_label='Média', line_cap="round", color='#e67e22', line_width=4)

p.grid.grid_line_alpha = 0
p.ygrid.band_fill_alpha = 0.1
p.ygrid.band_fill_color = "olive"
p.yaxis.axis_label = "Infectados"
p.xaxis.axis_label = "Geração"

show(p)